La cardiotocografia és una tècnica utilitzada per fer un monitoratge del ritme cardíac del fetus durant l’embaràs i el part. Es tracta d’una operació de monitoratge que, a part de la pulsació, pren informació de diverses variables per tal de detectar anomalies en el ritme del batec del cor, i determinar l’existència d’anomalies és per tant una tasca complicada de determinar de forma visual en temps real. Un dels aparells més utilitzats per fer aquest monitoratge és l’Omniview-SisPorto, un sistema automatitzat que realitza una anàlisi de la cardiotocografia i activa senyals d’alerta visual i sonora quan la combinació de mesures formen un patró considerat sospitós o patològic.
Malgrat el seu bon funcionament, s’han donat alguns casos en què patrons anòmals del ritme cardíac fetal han passat desapercebuts, i l’empresa Omniview ens ha contactat amb l’encàrrec de millorar els mètodes de classificació que fa el sistema SisPorto. L’objectiu és crear un model de classificació que millori els índexs actuals d’especificitat del sistema, de forma que es minimitzi el nombre de falsos negatius en la detecció d’anomalies del ritme cardíac fetal.
Per començar la nostra tasca, partim del dataset “Fetal cardiotocography data”, disponible a https://www.kaggle.com/datasets/akshat0007/fetalhr.
Es tracta d’un joc de dades amb 2129 observacions de 40 variables, que inclou dades capturades amb l’aparell SisPorto, però també altres variables recollides de forma manual per l’expert sanitari que realitzà la cardiotocografia en cada observació. Les dades inclouen també un diagnòstic que conclou, per part de l’expert mèdic, si hi ha un patró anòmal o potencialment patològic en cada observació. Aquest diagnòstic es reflecteix en dues de les variables del joc de dades. La primera, anomenada ‘Class’, és una escala de factors codificats de l’1 al 10, anant de normal a sospitós. La segona variable, també adient per aplicar algorismes de classificació, anomenada ‘NSP’, classifica les observacions en tres classes: 1=Normal, 2=Sospitós i 3=Patològic.
Presentem la totalitat de les variables que apareixen al joc de dades. En primer lloc tenim les variables que identifiquen les observacions.
-FileName (string): nom del fitxer de la cardiotocografia
-Data (string): data (mm/dd/aaaa) de l’exploració
-SegFile (string): identificador únic per a cada mesura.
Seguidament, tenim dues variables que determinen la duració (en segons) del monitoratge:
-b (numèrica): moment d’inici
-e (numèrica): moment de finalització
Les variables enregistrades pel sistema SisPorto:
-LB (numèrica): nivell (ritme cardíac) de base
-AC (numèrica): acceleracions
-FM (numèrica): moviments fetals
-UC (numèrica): contraccions uterines
-ASTV (numèrica): percentatge de temps amb variabilitat anormal de curta durada
-mSTV (numèrica): valor mitjà de variabilitat de curta durada
-ALTV (numèrica): percentatge de temps amb variabilitat anormal de llarga durada
-mLTV (numèrica): valor mitjà de variabilitat de llarga durada
Altres variables referents al ritme cardíac:
-LB (numèrica): nivell de base determinat per un expert
-DL (numèrica): desacceleracions lleugeres
-DS (numèrica): desacceleracions severes
-DP (numèrica): desacceleracions perllongades
-DR (numèrica): desacceleracions repetitives
Més enllà, el registre del ritme cardíac fetal genera un histograma, les dades del qual s’inclouen al dataset en les següents variables:
-Width (numèrica): amplada de l’histograma
-Min (numèrica): mínim de freqüència de l’histograma
-Max (numèrica): màxim de freqüència de l’histograma
-Nmax (numèrica): nombre de pics de l’histograma
-Nzeros (numèrica): nombre de zeros de l’histograma
-Mode (numèrica): mode
-Mean (numèrica): mitjana
-Median (numèrica): mediana
-Variance (numèrica): variància
-Tendency (categòrica): -1=asimètrica a l’esquerra, 0=simètrica, 1=asimètrica a la dreta
Finalment, hi ha un seguit de variables amb valoracions diagnòstiques de classificació:
-A (binària): son calmat
-B (binària): son REM
-C (binària): vetlla calmada
-D (binària): vetlla activa
-AD (binària): patró d’estrès
-DE (binària): patró de desacceleració
-LD (binària): patró sever de desacceleració
-FS (binària): patró sinusoïdal pla (estat patològic)
-SUSP (binària): patró sospitós
-CLASS (categòrica): codi de 0 a 1 per classes A a SUSP
-NSP (categòrica): 1=Normal, 2=Sospitós, 3=Patològic
Carreguem les dades:
CTG <- read.csv("CTG.csv")Mostrem les primeres línies del dataset:
as.data.frame(head(CTG))## FileName Date SegFile b e LBE LB AC FM UC ASTV MSTV ALTV
## 1 Variab10.txt 12/1/1996 CTG0001.txt 240 357 120 120 0 0 0 73 0.5 43
## 2 Fmcs_1.txt 5/3/1996 CTG0002.txt 5 632 132 132 4 0 4 17 2.1 0
## 3 Fmcs_1.txt 5/3/1996 CTG0003.txt 177 779 133 133 2 0 5 16 2.1 0
## 4 Fmcs_1.txt 5/3/1996 CTG0004.txt 411 1192 134 134 2 0 6 16 2.4 0
## 5 Fmcs_1.txt 5/3/1996 CTG0005.txt 533 1147 132 132 4 0 5 16 2.4 0
## 6 Fmcs_2.txt 5/3/1996 CTG0006.txt 0 953 134 134 1 0 10 26 5.9 0
## MLTV DL DS DP DR Width Min Max Nmax Nzeros Mode Mean Median Variance Tendency
## 1 2.4 0 0 0 0 64 62 126 2 0 120 137 121 73 1
## 2 10.4 2 0 0 0 130 68 198 6 1 141 136 140 12 0
## 3 13.4 2 0 0 0 130 68 198 5 1 141 135 138 13 0
## 4 23.0 2 0 0 0 117 53 170 11 0 137 134 137 13 1
## 5 19.9 0 0 0 0 117 53 170 9 0 137 136 138 11 1
## 6 0.0 9 0 2 0 150 50 200 5 3 76 107 107 170 0
## A B C D E AD DE LD FS SUSP CLASS NSP
## 1 0 0 0 0 0 0 0 0 1 0 9 2
## 2 0 0 0 0 0 1 0 0 0 0 6 1
## 3 0 0 0 0 0 1 0 0 0 0 6 1
## 4 0 0 0 0 0 1 0 0 0 0 6 1
## 5 0 1 0 0 0 0 0 0 0 0 2 1
## 6 0 0 0 0 0 0 0 1 0 0 8 3Verifiquem les dimensions del dataset:
dim(CTG)## [1] 2129 40
Mostrem un resum de les variables:
summary(CTG)## FileName Date SegFile b
## Length:2129 Length:2129 Length:2129 Min. : 0.0
## Class :character Class :character Class :character 1st Qu.: 55.0
## Mode :character Mode :character Mode :character Median : 538.0
## Mean : 878.4
## 3rd Qu.:1521.0
## Max. :3296.0
## NA's :3
## e LBE LB AC
## Min. : 287 Min. :106.0 Min. :106.0 Min. : 0.000
## 1st Qu.:1009 1st Qu.:126.0 1st Qu.:126.0 1st Qu.: 0.000
## Median :1241 Median :133.0 Median :133.0 Median : 1.000
## Mean :1703 Mean :133.3 Mean :133.3 Mean : 2.722
## 3rd Qu.:2435 3rd Qu.:140.0 3rd Qu.:140.0 3rd Qu.: 4.000
## Max. :3599 Max. :160.0 Max. :160.0 Max. :26.000
## NA's :3 NA's :3 NA's :3 NA's :3
## FM UC ASTV MSTV
## Min. : 0.000 Min. : 0.000 Min. :12.00 Min. :0.200
## 1st Qu.: 0.000 1st Qu.: 1.000 1st Qu.:32.00 1st Qu.:0.700
## Median : 0.000 Median : 3.000 Median :49.00 Median :1.200
## Mean : 7.503 Mean : 3.669 Mean :47.01 Mean :1.335
## 3rd Qu.: 2.000 3rd Qu.: 5.000 3rd Qu.:61.00 3rd Qu.:1.700
## Max. :564.000 Max. :23.000 Max. :87.00 Max. :7.000
## NA's :2 NA's :2 NA's :2 NA's :2
## ALTV MLTV DL DS
## Min. : 0.000 Min. : 0.000 Min. : 0.000 Min. :0.000000
## 1st Qu.: 0.000 1st Qu.: 4.600 1st Qu.: 0.000 1st Qu.:0.000000
## Median : 0.000 Median : 7.400 Median : 0.000 Median :0.000000
## Mean : 9.885 Mean : 8.208 Mean : 1.576 Mean :0.003759
## 3rd Qu.:11.000 3rd Qu.:10.800 3rd Qu.: 3.000 3rd Qu.:0.000000
## Max. :91.000 Max. :50.700 Max. :16.000 Max. :1.000000
## NA's :2 NA's :2 NA's :1 NA's :1
## DP DR Width Min Max
## Min. :0.0000 Min. :0 Min. : 3.00 Min. : 50.00 Min. :122
## 1st Qu.:0.0000 1st Qu.:0 1st Qu.: 37.00 1st Qu.: 67.00 1st Qu.:152
## Median :0.0000 Median :0 Median : 67.50 Median : 93.00 Median :162
## Mean :0.1278 Mean :0 Mean : 70.45 Mean : 93.58 Mean :164
## 3rd Qu.:0.0000 3rd Qu.:0 3rd Qu.:100.00 3rd Qu.:120.00 3rd Qu.:174
## Max. :4.0000 Max. :0 Max. :180.00 Max. :159.00 Max. :238
## NA's :1 NA's :1 NA's :3 NA's :3 NA's :3
## Nmax Nzeros Mode Mean
## Min. : 0.000 Min. : 0.0000 Min. : 60.0 Min. : 73.0
## 1st Qu.: 2.000 1st Qu.: 0.0000 1st Qu.:129.0 1st Qu.:125.0
## Median : 3.000 Median : 0.0000 Median :139.0 Median :136.0
## Mean : 4.068 Mean : 0.3236 Mean :137.5 Mean :134.6
## 3rd Qu.: 6.000 3rd Qu.: 0.0000 3rd Qu.:148.0 3rd Qu.:145.0
## Max. :18.000 Max. :10.0000 Max. :187.0 Max. :182.0
## NA's :3 NA's :3 NA's :3 NA's :3
## Median Variance Tendency A
## Min. : 77.0 Min. : 0.00 Min. :-1.0000 Min. :0.0000
## 1st Qu.:129.0 1st Qu.: 2.00 1st Qu.: 0.0000 1st Qu.:0.0000
## Median :139.0 Median : 7.00 Median : 0.0000 Median :0.0000
## Mean :138.1 Mean : 18.81 Mean : 0.3203 Mean :0.1806
## 3rd Qu.:148.0 3rd Qu.: 24.00 3rd Qu.: 1.0000 3rd Qu.:0.0000
## Max. :186.0 Max. :269.00 Max. : 1.0000 Max. :1.0000
## NA's :3 NA's :3 NA's :3 NA's :3
## B C D E
## Min. :0.0000 Min. :0.00000 Min. :0.0000 Min. :0.00000
## 1st Qu.:0.0000 1st Qu.:0.00000 1st Qu.:0.0000 1st Qu.:0.00000
## Median :0.0000 Median :0.00000 Median :0.0000 Median :0.00000
## Mean :0.2723 Mean :0.02493 Mean :0.0381 Mean :0.03387
## 3rd Qu.:1.0000 3rd Qu.:0.00000 3rd Qu.:0.0000 3rd Qu.:0.00000
## Max. :1.0000 Max. :1.00000 Max. :1.0000 Max. :1.00000
## NA's :3 NA's :3 NA's :3 NA's :3
## AD DE LD FS
## Min. :0.0000 Min. :0.0000 Min. :0.00000 Min. :0.00000
## 1st Qu.:0.0000 1st Qu.:0.0000 1st Qu.:0.00000 1st Qu.:0.00000
## Median :0.0000 Median :0.0000 Median :0.00000 Median :0.00000
## Mean :0.1562 Mean :0.1185 Mean :0.05033 Mean :0.03246
## 3rd Qu.:0.0000 3rd Qu.:0.0000 3rd Qu.:0.00000 3rd Qu.:0.00000
## Max. :1.0000 Max. :1.0000 Max. :1.00000 Max. :1.00000
## NA's :3 NA's :3 NA's :3 NA's :3
## SUSP CLASS NSP
## Min. :0.00000 Min. : 1.00 Min. :1.000
## 1st Qu.:0.00000 1st Qu.: 2.00 1st Qu.:1.000
## Median :0.00000 Median : 4.00 Median :1.000
## Mean :0.09266 Mean : 4.51 Mean :1.304
## 3rd Qu.:0.00000 3rd Qu.: 7.00 3rd Qu.:1.000
## Max. :1.00000 Max. :10.00 Max. :3.000
## NA's :3 NA's :3 NA's :3
Si seleccionem les files amb valors NA veiem que es tracta de les últimes tres files del dataset, on la majoria de camps estan buits:
CTG[rowSums(is.na(CTG)) > 0,]## FileName Date SegFile b e LBE LB AC FM UC ASTV MSTV ALTV MLTV DL DS DP
## 2127 NA NA NA NA NA NA NA NA NA NA NA NA NA NA
## 2128 NA NA NA NA NA NA NA NA NA NA NA 0 0 0
## 2129 NA NA NA NA NA 564 23 87 7 91 50.7 16 1 4
## DR Width Min Max Nmax Nzeros Mode Mean Median Variance Tendency A B C
## 2127 NA NA NA NA NA NA NA NA NA NA NA NA NA NA
## 2128 0 NA NA NA NA NA NA NA NA NA NA NA NA NA
## 2129 0 NA NA NA NA NA NA NA NA NA NA NA NA NA
## D E AD DE LD FS SUSP CLASS NSP
## 2127 NA NA NA NA NA NA NA NA NA
## 2128 NA NA NA NA NA NA NA NA NA
## 2129 NA NA NA NA NA NA NA NA NA
Així doncs, eliminem aquests darrers tres registres sense perdre informació valuosa:
CTG2 <- na.omit(CTG)
dim(CTG2)## [1] 2126 40
Comprovem que amb aquesta operació, el dataset queda net de valors buits:
sum(is.na(CTG2))## [1] 0Seguidament, observem que hi ha dos variables que no ens semblen d’especial rellevància de cara a la tasca. Es tracta de les variables ‘FileName’ i ‘Date’. Considerem, d’una banda, que els registres queden correctament identificats si mantenim la variable ‘SegFile’, i d’altra banda, la data del registre no hauria de tenir cap influència en la classificació (assumint que la precisió de l’instrument de mesura no es veu afectat per la data de la mesura). Així doncs, eliminem aquestes dues variables del dataset:
CTG2 <- subset(CTG2, select = -c(FileName, Date))Seguidament, obtenint informació suplementària sobre el joc de dades al repositori de UCI, ens adonem que les variables ‘AC’, ‘FM’, ‘UC’, ‘DL’, ‘DS’ i ‘DP’ apareixen al dataset com a nombres absoluts d’ocurrències. Per tal que aquestes variables es puguin comparar entre registres, cal que estiguin expressades en nombre d’ocurrències per segon. Com que cada registre té una durada diferent, cal fer aquesta operació registre a registre. Per fer-ho, utilitzem les variables ‘b’ i ‘e’ (moment d’inici i de finalització) per determinar la durada del registre, i creem noves variables que expressen ocurrències per segon:
CTG2$ACps <- round(CTG2$AC / (CTG2$e - CTG2$b), digits = 3)
CTG2$FMps <- round(CTG2$FM / (CTG2$e - CTG2$b), digits = 3)
CTG2$UCps <- round(CTG2$UC / (CTG2$e - CTG2$b), digits = 3)
CTG2$DLps <- round(CTG2$DL / (CTG2$e - CTG2$b), digits = 3)
CTG2$DSps <- round(CTG2$DS / (CTG2$e - CTG2$b), digits = 3)
CTG2$DPps <- round(CTG2$DP / (CTG2$e - CTG2$b), digits = 3)Mostrem un nou resum per seguir explorant les dades:
str(CTG2)## 'data.frame': 2126 obs. of 44 variables:
## $ SegFile : chr "CTG0001.txt" "CTG0002.txt" "CTG0003.txt" "CTG0004.txt" ...
## $ b : int 240 5 177 411 533 0 240 62 120 181 ...
## $ e : int 357 632 779 1192 1147 953 953 679 779 1192 ...
## $ LBE : int 120 132 133 134 132 134 134 122 122 122 ...
## $ LB : int 120 132 133 134 132 134 134 122 122 122 ...
## $ AC : int 0 4 2 2 4 1 1 0 0 0 ...
## $ FM : int 0 0 0 0 0 0 0 0 0 0 ...
## $ UC : int 0 4 5 6 5 10 9 0 1 3 ...
## $ ASTV : int 73 17 16 16 16 26 29 83 84 86 ...
## $ MSTV : num 0.5 2.1 2.1 2.4 2.4 5.9 6.3 0.5 0.5 0.3 ...
## $ ALTV : int 43 0 0 0 0 0 0 6 5 6 ...
## $ MLTV : num 2.4 10.4 13.4 23 19.9 0 0 15.6 13.6 10.6 ...
## $ DL : int 0 2 2 2 0 9 6 0 0 0 ...
## $ DS : int 0 0 0 0 0 0 0 0 0 0 ...
## $ DP : int 0 0 0 0 0 2 2 0 0 0 ...
## $ DR : int 0 0 0 0 0 0 0 0 0 0 ...
## $ Width : int 64 130 130 117 117 150 150 68 68 68 ...
## $ Min : int 62 68 68 53 53 50 50 62 62 62 ...
## $ Max : int 126 198 198 170 170 200 200 130 130 130 ...
## $ Nmax : int 2 6 5 11 9 5 6 0 0 1 ...
## $ Nzeros : int 0 1 1 0 0 3 3 0 0 0 ...
## $ Mode : int 120 141 141 137 137 76 71 122 122 122 ...
## $ Mean : int 137 136 135 134 136 107 107 122 122 122 ...
## $ Median : int 121 140 138 137 138 107 106 123 123 123 ...
## $ Variance: int 73 12 13 13 11 170 215 3 3 1 ...
## $ Tendency: int 1 0 0 1 1 0 0 1 1 1 ...
## $ A : int 0 0 0 0 0 0 0 0 0 0 ...
## $ B : int 0 0 0 0 1 0 0 0 0 0 ...
## $ C : int 0 0 0 0 0 0 0 0 0 0 ...
## $ D : int 0 0 0 0 0 0 0 0 0 0 ...
## $ E : int 0 0 0 0 0 0 0 0 0 0 ...
## $ AD : int 0 1 1 1 0 0 0 0 0 0 ...
## $ DE : int 0 0 0 0 0 0 0 0 0 0 ...
## $ LD : int 0 0 0 0 0 1 1 0 0 0 ...
## $ FS : int 1 0 0 0 0 0 0 1 1 1 ...
## $ SUSP : int 0 0 0 0 0 0 0 0 0 0 ...
## $ CLASS : int 9 6 6 6 2 8 8 9 9 9 ...
## $ NSP : int 2 1 1 1 1 3 3 3 3 3 ...
## $ ACps : num 0 0.006 0.003 0.003 0.007 0.001 0.001 0 0 0 ...
## $ FMps : num 0 0 0 0 0 0 0 0 0 0 ...
## $ UCps : num 0 0.006 0.008 0.008 0.008 0.01 0.013 0 0.002 0.003 ...
## $ DLps : num 0 0.003 0.003 0.003 0 0.009 0.008 0 0 0 ...
## $ DSps : num 0 0 0 0 0 0 0 0 0 0 ...
## $ DPps : num 0 0 0 0 0 0.002 0.003 0 0 0 ...Tenim la impressió que pot ser que no hi hagi diferència entre la línia de base mesurada pel sistema SisPorto i la línia de base establerta per l’expert:
all(CTG2$LBE == CTG2$LB)## [1] TRUEEfectivament, les dues variables són idèntiques, i per tant en podem eliminar una de les dues:
CTG2 <- subset(CTG2, select = -LBE)Ens adonem també que tots els valors de la variable ‘DR’ són 0, i per tant també eliminem aquesta variable:
summary(CTG2$DR)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0 0 0 0 0 0CTG3 <- subset(CTG2, select = -DR)Hem analitzat una per una les variables numèriques que considerem interessants de cara a una possible tasca de classificació. A continuació mostrem algunes de les variables on hem detectat valors aïllats:
La variable ‘AC’ (acceleracions), mostra un gran nombre de valors aïllats:
boxplot.stats(CTG3$AC)$out## [1] 12 11 13 12 11 15 11 13 12 17 13 11 17 11 18 12 11 11 13 16 14 14 12 11 13
## [26] 12 14 12 14 12 12 11 12 12 19 19 13 17 13 13 15 12 14 26 18 12 17 13 14 13
## [51] 17 14 14 14 11 11 12 13 11 11 14 13 14 11 15 13 15 14 13 12 15 14 16 13 12
## [76] 17 17 16 16 11 12 11 21boxplot(CTG3$AC, main = "AC")
Si mostrem el valor més allunyat, veiem que la resta de variables no mostren valors estranys, i per tant la totalitat del registre sembla tenir coherència. Es podria tractar d’un error només en el valor de la variable ‘AC’. Fixant-nos-hi més en detall, veiem que la duració del registre és de més de 3000 segons, una duració molt llarga, que fa que un nombre total gran d’ocurrències de la variable ‘AC’ no sigui tan estrany:
CTG3[CTG3$AC > 25,]## SegFile b e LB AC FM UC ASTV MSTV ALTV MLTV DL DS DP Width Min
## 1648 CTG1650.txt 357 3591 130 26 3 13 51 1.6 0 4.2 12 0 0 109 77
## Max Nmax Nzeros Mode Mean Median Variance Tendency A B C D E AD DE LD FS
## 1648 186 4 0 144 138 143 38 0 0 0 0 0 0 1 0 0 0
## SUSP CLASS NSP ACps FMps UCps DLps DSps DPps
## 1648 0 6 1 0.008 0.001 0.004 0.004 0 0Per comparar-ho, fem la mateixa operació de detecció d’outliers, però aquest cop utilitzem la variable ‘ACps’, que mostra ocurrències per segon en lloc de nombre absolut d’ocurrències:
boxplot.stats(CTG3$ACps)$out## [1] 0.017 0.016 0.019 0.016 0.016 0.016 0.017 0.016 0.018 0.017 0.018 0.016
## [13] 0.017 0.016boxplot(CTG3$ACps, main = "ACps")
Veiem que els valors extrems són menys, i també menys allunyats. Considerem que aquests valors allunyats no són errors de registre, i els mantenim en el dataset sense cap correcció. Realitzem la resta de l’anàlisi de detecció d’outliers utilitzant les variables que registren ocurrències per segon:
CTG3[CTG3$ACps > 0.016,]## SegFile b e LB AC FM UC ASTV MSTV ALTV MLTV DL DS DP Width Min
## 182 CTG0182.txt 194 948 138 13 0 4 35 5.3 0 4.5 0 0 0 148 52
## 530 CTG0530.txt 829 1192 142 7 31 0 32 2.3 0 0.0 0 0 0 144 56
## 631 CTG0631.txt 136 941 134 14 2 3 48 2.2 0 0.0 0 0 0 120 50
## 1095 CTG1095.txt 1548 2114 122 10 0 1 22 2.5 0 2.2 0 0 0 52 100
## 1097 CTG1097.txt 1352 1871 123 9 0 1 24 2.2 0 1.7 0 0 0 52 100
## 1249 CTG1250.txt 2337 2897 112 10 0 4 25 1.4 0 0.6 0 0 0 35 109
## 1860 CTG1862.txt 828 1648 138 14 0 3 51 0.9 0 0.2 0 0 0 49 122
## Max Nmax Nzeros Mode Mean Median Variance Tendency A B C D E AD DE LD FS
## 182 200 11 2 146 157 161 72 1 0 0 0 1 0 0 0 0 0
## 530 200 10 0 170 158 162 37 1 0 0 0 1 0 0 0 0 0
## 631 170 5 0 160 150 155 28 1 0 0 0 1 0 0 0 0 0
## 1095 152 1 0 136 132 134 6 0 0 1 0 0 0 0 0 0 0
## 1097 152 1 0 136 133 135 5 0 0 1 0 0 0 0 0 0 0
## 1249 144 0 0 116 120 119 6 -1 0 1 0 0 0 0 0 0 0
## 1860 171 4 0 147 149 150 5 0 0 1 0 0 0 0 0 0 0
## SUSP CLASS NSP ACps FMps UCps DLps DSps DPps
## 182 0 4 1 0.017 0.000 0.005 0 0 0
## 530 0 4 1 0.019 0.085 0.000 0 0 0
## 631 0 4 1 0.017 0.002 0.004 0 0 0
## 1095 0 2 1 0.018 0.000 0.002 0 0 0
## 1097 0 2 1 0.017 0.000 0.002 0 0 0
## 1249 0 2 1 0.018 0.000 0.007 0 0 0
## 1860 0 2 1 0.017 0.000 0.004 0 0 0Analitzem la variable ‘FMps’ (moviments fetals per segon):
boxplot.stats(CTG3$FMps)$out## [1] 0.072 0.222 0.408 0.380 0.441 0.383 0.451 0.469 0.340 0.425 0.334 0.135
## [13] 0.099 0.108 0.112 0.089 0.103 0.085 0.109 0.079 0.065 0.055 0.058 0.047
## [25] 0.038 0.012 0.018 0.020 0.009 0.010 0.008 0.008 0.028 0.026 0.107 0.009
## [37] 0.008 0.013 0.016 0.029 0.050 0.053 0.058 0.009 0.011 0.008 0.013 0.013
## [49] 0.016 0.008 0.015 0.008 0.010 0.022 0.026 0.015 0.015 0.010 0.010 0.021
## [61] 0.008 0.016 0.018 0.013 0.017 0.016 0.019 0.008 0.009 0.009 0.010 0.009
## [73] 0.008 0.009 0.010 0.009 0.013 0.013 0.018 0.015 0.012 0.019 0.025 0.012
## [85] 0.020 0.018 0.020 0.021 0.020 0.014 0.021 0.022 0.025 0.022 0.024 0.028
## [97] 0.016 0.015 0.019 0.024 0.009 0.019 0.011 0.016 0.019 0.021 0.010 0.019
## [109] 0.008 0.019 0.022 0.009 0.009 0.009 0.009 0.009 0.011 0.013 0.013 0.016
## [121] 0.019 0.013 0.008 0.008 0.009 0.012 0.009 0.017 0.023 0.023 0.020 0.020
## [133] 0.024 0.035 0.054 0.030 0.048 0.088 0.030 0.054 0.043 0.009 0.054 0.011
## [145] 0.058 0.052 0.085 0.091 0.026 0.033 0.092 0.084 0.115 0.084 0.015 0.017
## [157] 0.013 0.014 0.014 0.015 0.035 0.033 0.029 0.041 0.024 0.040 0.027 0.029
## [169] 0.031 0.026 0.063 0.060 0.050 0.085 0.071 0.060 0.029 0.008 0.008 0.008
## [181] 0.020 0.010 0.027 0.017 0.018 0.010 0.008 0.010 0.009 0.014 0.010 0.013
## [193] 0.025 0.009 0.013 0.017 0.017 0.021 0.025 0.018 0.013 0.017 0.010 0.008
## [205] 0.306 0.298 0.139 0.189 0.014 0.157 0.235 0.408 0.360 0.455 0.443 0.470
## [217] 0.477 0.446 0.481 0.369 0.335 0.430 0.346 0.323 0.375 0.346 0.353 0.010
## [229] 0.021 0.017 0.045 0.016 0.011 0.008 0.017 0.017 0.020 0.030 0.012 0.043
## [241] 0.035 0.040 0.048 0.012 0.014 0.013 0.011 0.008 0.010 0.010 0.012 0.012
## [253] 0.012 0.014 0.008 0.010 0.012 0.009 0.011 0.010 0.014 0.011 0.010 0.010
## [265] 0.010 0.010 0.010 0.011 0.010 0.016 0.010 0.008 0.009 0.008 0.011 0.011
## [277] 0.017 0.019 0.019 0.025 0.028 0.048 0.043 0.032 0.040 0.051 0.036 0.037
## [289] 0.050 0.035 0.049 0.049 0.041 0.050 0.048 0.054 0.052 0.008 0.008 0.015
## [301] 0.017 0.008 0.009 0.009 0.009 0.010 0.009boxplot(CTG3$FMps, main = "FMps")
hist(CTG3$FMps, breaks = sqrt(nrow(CTG3)))
En el cas de la variable FMps, hi ha en un principi una quantitat extrema d’outliers. L’anàlisi visual mostra que en realitat no es tracta d’outliers, sinó que estem davant d’una variable amb una distribució molt asimètrica. La variable ‘FM’ mesura ‘moviments fetals’, i sembla lògic que hi hagi una gran variabilitat, així com una gran quantitat de valors baixos o zero, ja que quan el fetus dorm o està en estat de vetlla calmada, els moviments són pocs. En tot cas, creem una nova variable amb una transformació logarítmica d’aquesta variable, en cas que ens pugui resultar útil en models posteriors:
library(ggplot2)
CTG3$FMps_log <- log(CTG3$FMps)
ggplot(mapping = aes(CTG3$FMps_log)) + geom_density()## Warning: Removed 1311 rows containing non-finite values (stat_density).
Analitzem la variable ‘UCps’(contraccions uterines per segon):
boxplot.stats(CTG3$UCps)$out## [1] 0.015boxplot(CTG3$UCps, main = "UCps")
Detectem un únic valor allunyat, però no tenim evidència per considerar-lo com a erroni, i per tant el mantenim en el dataset:
CTG3[CTG3$UCps > 0.014,]## SegFile b e LB AC FM UC ASTV MSTV ALTV MLTV DL DS DP Width Min
## 1165 CTG1165.txt 910 1379 131 5 0 7 26 1.5 0 3 0 0 0 61 109
## Max Nmax Nzeros Mode Mean Median Variance Tendency A B C D E AD DE LD FS
## 1165 170 2 1 155 151 154 11 1 0 0 0 1 0 0 0 0 0
## SUSP CLASS NSP ACps FMps UCps DLps DSps DPps FMps_log
## 1165 0 4 1 0.011 0 0.015 0 0 0 -Inf
Analitzem la variable ‘ASTV’ (percentatge de temps amb variabilitat anormal de curta durada), i no hi trobem outliers:
boxplot.stats(CTG3$ASTV)$out## integer(0)
Analitzem la variable ‘ALTV’ (percentatge de temps amb variabilitat anormal de llarga durada), i ens trobem amb la mateixa situació que amb la variable ‘FM’, una distribució molt asimètrica, amb molts valors allunyats, però que no es poden considerar com a registres erronis. No hi fem cap modificació:
boxplot.stats(CTG3$ALTV)$out## [1] 43 79 72 71 40 69 54 53 38 29 37 29 67 68 75 74 30 49 39 31 32 34 38 31 58
## [26] 39 46 57 75 29 33 40 51 37 34 34 54 58 52 32 38 59 78 84 71 59 52 62 39 45
## [51] 30 45 35 67 72 45 58 32 41 61 56 84 56 84 44 39 71 71 61 67 81 42 57 56 61
## [76] 77 67 88 91 78 84 47 41 59 55 61 49 49 82 86 29 67 78 84 60 68 30 32 41 28
## [101] 38 29 49 35 62 68 70 47 48 51 31 43 50 47 71 48 39 39 38 31 31 34 28 29 38
## [126] 32 37 32 37 36 28 28 31 40 72 85 64 58 66 91 90 54 70 77 58 81 29 40 38 54
## [151] 73 43 50 42 53 64 41 34 33 34 37 44 62 29 69 72 91 90 84 42 91 66 40 34 32
## [176] 30 32 38 32 45 32 32 32 31 33 39 39 44 55 41 33 31 32 32 35 53 44 31 33 57
## [201] 50 64 59 67 62 71 74 59 48 45 62 59 50 61 39 43 56 62 53 44 48 61 54 37 58
## [226] 40 30 37 42 28 47 48 28 43 38 28 39 32 61 48 35 34 45 49 42 33 32 32 31 33
## [251] 35 32 31 28 33 36 49 52 46 41 37 54 60 62 59 62 56 49 38 32 45 30 44 30 46
## [276] 37 40 54 65 46 34 44 46 55 58 59 64 60 61 73 71 67 63 52 58 35 28 73 70 59
## [301] 44 34 29 40 36 33 48 36 36boxplot(CTG3$ALTV, main = "ALTV")
Pel que fa a la variable ‘DS’, ens trobem amb una situació diferent. Els possibles valors són 0 i 1, i dels 2126 casos, només 7 es classifiquen com a ‘1’. La variable mesura ‘desacceleracions severes’, i és per tant normal que només una petita minoria dels casos presentin aquest patró:
CTG3[CTG3$DS != 0,]## SegFile b e LB AC FM UC ASTV MSTV ALTV MLTV DL DS DP Width Min
## 1489 CTG1491.txt 1532 2655 132 2 0 9 31 1.4 0 11.5 0 1 1 102 61
## 1490 CTG1492.txt 1728 2655 132 0 0 6 32 1.3 0 13.6 0 1 1 91 60
## 1792 CTG1794.txt 1839 2894 121 0 1 4 66 2.1 0 6.4 11 1 0 105 55
## 1793 CTG1795.txt 1922 2894 121 0 1 3 67 2.1 0 0.0 11 1 0 102 55
## 1794 CTG1796.txt 1922 2771 121 0 1 4 66 2.1 0 0.0 10 1 0 102 55
## 1795 CTG1797.txt 2018 2892 121 0 1 3 68 2.1 0 0.0 9 1 0 102 55
## 1796 CTG1798.txt 2153 2892 121 0 0 3 70 1.9 0 0.0 7 1 0 102 55
## Max Nmax Nzeros Mode Mean Median Variance Tendency A B C D E AD DE LD FS
## 1489 163 5 0 99 121 129 94 1 0 0 0 0 0 1 0 0 0
## 1490 151 1 1 99 116 125 72 1 0 0 0 0 0 0 0 1 0
## 1792 160 7 0 67 85 92 109 -1 0 0 0 0 0 0 0 1 0
## 1793 157 4 1 67 81 87 89 -1 0 0 0 0 0 0 0 1 0
## 1794 157 5 1 67 83 90 98 -1 0 0 0 0 0 0 0 1 0
## 1795 157 3 1 67 79 82 83 -1 0 0 0 0 0 0 0 1 0
## 1796 157 6 2 67 76 79 68 -1 0 0 0 0 0 0 0 1 0
## SUSP CLASS NSP ACps FMps UCps DLps DSps DPps FMps_log
## 1489 0 6 1 0.002 0.000 0.008 0.000 0.001 0.001 -Inf
## 1490 0 8 3 0.000 0.000 0.006 0.000 0.001 0.001 -Inf
## 1792 0 8 3 0.000 0.001 0.004 0.010 0.001 0.000 -6.907755
## 1793 0 8 3 0.000 0.001 0.003 0.011 0.001 0.000 -6.907755
## 1794 0 8 3 0.000 0.001 0.005 0.012 0.001 0.000 -6.907755
## 1795 0 8 3 0.000 0.001 0.003 0.010 0.001 0.000 -6.907755
## 1796 0 8 3 0.000 0.000 0.004 0.009 0.001 0.000 -InfProcedim a discretitzar diverses de les variables. En primer lloc, afegirem etiquetes a la variable ‘DS’ que acabem d’analitzar, creant una nova variable amb els codis ‘absent’ i ‘present’:
library(arules)## Loading required package: Matrix##
## Attaching package: 'arules'## The following objects are masked from 'package:base':
##
## abbreviate, writeCTG3$DS_d <- discretize(CTG3$DS, method = "fixed", breaks = c(-Inf, 1, Inf),labels = c("absent", "present"))
table(CTG3$DS_d)##
## absent present
## 2119 7Seguidament discretitzarem totes les variables que mesuren ocurrències per segon, creant en cada cas tres nivells: ‘baix’, ‘mitjà’ i ‘alt’. Farem la discretització utilitzant el mètode k-means:
CTG3$AC_d <- discretize(CTG3$ACps, method = "cluster", labels = c("low", "medium", "high"))
CTG3$FM_d <- discretize(CTG3$FMps, method = "cluster", labels = c("low", "medium", "high"))
CTG3$UC_d <- discretize(CTG3$UCps, method = "cluster", labels = c("low", "medium", "high"))
CTG3$DL_d <- discretize(CTG3$DLps, method = "cluster", labels = c("low", "medium", "high"))
CTG3$DP_d <- discretize(CTG3$DPps, method = "cluster", labels = c("low", "medium", "high"))
summary(CTG3$AC_d)## low medium high
## 1037 660 429summary(CTG3$FM_d)## low medium high
## 2037 61 28summary(CTG3$UC_d)## low medium high
## 450 616 1060summary(CTG3$DL_d)## low medium high
## 1627 404 95summary(CTG3$DP_d)## low medium high
## 1948 70 108
També crearem una nova variable ‘NSP2’ a partir de la variable ‘NSP’, on ajuntarem les categories ‘sospitós’ i ‘patològic’. Amb això obtenim una variable target classificatòria binària, i podem més endavant triar si els nostres models de classificació tenen 2 o 3 classes a identificar:
summary(as.factor(CTG3$NSP))## 1 2 3
## 1655 295 176CTG3$NSP2 <- as.factor(CTG3$NSP)
levels(CTG3$NSP2) <- c('Normal', 'Suspect', 'Pathologic')
levels(CTG3$NSP2)[levels(CTG3$NSP2)=="Pathologic"] <-"Suspect"
summary(CTG3$NSP2)## Normal Suspect
## 1655 471Finalment, i abans de donar per acabat el treball de neteja de les dades, ens assegurem que totes les variables categòriques o binàries estiguin factoritzades:
factcols <- c('A', 'B', 'C', 'D', 'AD', 'DE', 'LD', 'FS', 'SUSP', 'CLASS', 'NSP')
CTG3[factcols] <- lapply(CTG3[factcols], factor)
sapply(CTG3, class)## SegFile b e LB AC FM
## "character" "integer" "integer" "integer" "integer" "integer"
## UC ASTV MSTV ALTV MLTV DL
## "integer" "integer" "numeric" "integer" "numeric" "integer"
## DS DP Width Min Max Nmax
## "integer" "integer" "integer" "integer" "integer" "integer"
## Nzeros Mode Mean Median Variance Tendency
## "integer" "integer" "integer" "integer" "integer" "integer"
## A B C D E AD
## "factor" "factor" "factor" "factor" "integer" "factor"
## DE LD FS SUSP CLASS NSP
## "factor" "factor" "factor" "factor" "factor" "factor"
## ACps FMps UCps DLps DSps DPps
## "numeric" "numeric" "numeric" "numeric" "numeric" "numeric"
## FMps_log DS_d AC_d FM_d UC_d DL_d
## "numeric" "factor" "factor" "factor" "factor" "factor"
## DP_d NSP2
## "factor" "factor"Creem un nou dataframe amb totes les transformacions realitzades:
ctg1 <- CTG3
write.csv(ctg1, "ctg1.csv", row.names = FALSE)ctgp2 <- read.csv("ctg1.csv")Com vam comentar a la pràctica anterior, l’objectiu del nostre estudi és millorar la capacitat de classificació de models construïts amb les dades proporcionades per sistema SisPorto. El nostre joc de dades conté una gran quantitat de dades, tant numèriques com categòriques, però seguint el nostre objectiu, seleccionem les variables numèriques que genera els sistema SisPorto per aplicar-hi un model no supervisat. Prescindim de DSps (desacceleracions severes), ja que entenem que quan es detecta una desacceleració severa, això per si sol ja indica un cas patològic. Ens interessa en canvi veure com les altres variables, que capturen anomalies més lleus o subtils, es combinen dintre de cada grup diagnòstic.
num_df <- ctgp2[c('ACps', 'FMps', 'UCps', 'DLps', 'DPps', 'ASTV', 'ALTV')]Normalitzem les dades:
#definim la funció de normalització
mm_norm <- function(x) {
(x - min(x)) / (max(x) - min(x))
}
num_norm <- as.data.frame(lapply(num_df, mm_norm))
if (!require('cluster')) install.packages('cluster')## Loading required package: clusterlibrary(cluster)Per tal de començar a trobar el nombre ideal de clústers en un model no supervisat, iterem la funció kmeans per agrupacions d’entre 2 i 10 clústers:
d <- daisy(num_norm)
results <- rep(0, 10)
for (i in c(2,3,4,5,6,7,8,9,10))
{
fit <- kmeans(num_norm, i)
y_cluster <- fit$cluster
sk <- silhouette(y_cluster, d)
results[i] <- mean(sk[,3])
}
Prenent com a mesura l’índex de silhouette, el nombre òptim de clústers seria 3, que precisament coincideix amb els tres grups que tenim definits segons la variable NSP (‘normal’, ‘sospitós’ i ‘patològic’):
plot(2:10,results[2:10],type="o",col="blue",pch=0,xlab="Number of clusters",ylab="Silhouette")
Si, com a alternativa, fem el càlcul de la suma de quadrats de distàncies respecte a centroides, també per a entre 2 i 10 agrupacions, el gràfic resultant és més difícil d’interpretar que l’anterior, però sembla que hi pugui haver un ‘colze’ també a l’alçada de 3 clústers:
results <- rep(0, 10)
for (i in c(2,3,4,5,6,7,8,9,10))
{
fit <- kmeans(num_norm, i)
results[i] <- fit$tot.withinss
}
plot(2:10,results[2:10],type="o",col="blue",pch=0,xlab="Nombre de clústers",ylab="tot.tot.withinss")
Pel que fa als índexs ‘Calinski-Harabasz’ (‘ch’) i silueta mitjana (‘asw’), aquests ens dónen respectivament un nombre òptim de clústers de 2 i 3:
if (!require('fpc')) install.packages('fpc'); library('fpc')## Loading required package: fpcfit_ch <- kmeansruns(num_norm, krange = 1:10, criterion = "ch")
fit_asw <- kmeansruns(num_norm, krange = 1:10, criterion = "asw")
fit_ch$bestk## [1] 2fit_asw$bestk## [1] 3A partir d’aquest moment, per avaluar l’algorisme de classificació, tenim dues opcions. Podem comparar la classificació amb 3 clústers comparant-la amb la variable NSP, que determina 3 nivells, o també podem crear un model amb 2 clústers, i comparar-lo amb la variable NSP2 que havíem creat anteriorment, i que ajuntava els nivells ‘sospitós’ i ‘patològic’ en un de sol.
Veiem una primera comparació entre les variables ALTV i ASTV amb 3 clústers:
ctg3clusters <- kmeans(num_norm, 3)
plot(num_norm[c(6,7)], col=ctg3clusters$cluster, main="Classificació k-means")
Com podem veure, si comparem aquest gràfic amb la classificació real, els tres clústers creats per l’algorisme no corresponen de forma gaire fidel als 3 grups diagnòstics. Hi ha una mena d’inversió dels grups. El grup definit pels nivells baixos d’ambdues variables en la classificació real (en negre al gràfic de sota) es veu dividit en 2 subgrups en la classificació amb k-means, i els dos altres grups reals, en la zona de valors progressivament més alts, es veuen compactats en un sol grup.
plot(num_norm[c(6,7)], col=as.factor(ctgp2$NSP), main="Classificació real")
Provem ara de comparar la classificació de k-means de 2 clústers amb la variable NSP2, que utilitza només 2 etiquetes diagnòstiques.
ctg2clusters <- kmeans(num_norm, 2)
plot(num_norm[c(6,7)], col=ctg2clusters$cluster, main="Classificació k-means")
plot(num_norm[c(6,7)], col=as.factor(ctgp2$NSP2), main="Classificació real")
A la comparació dels gràfics veiem que, a grans trets, la classificació sembla coincidir, si més no pels valors extrems de les dues variables ALTV i ASTV. Però si ens hi fixem, en la zona intermèdia, la pertinença a un grup o l’altre difereix notablement entre les dues gràfiques. Mentre que a la classificació real l’encavalcament és més gradual, i fins i tot hi ha casos pertanyents al primer grup que s’enfilen fins els valors més alts de les dues variables, a la classificació amb k-means els dos grups es veuen dividits sense pràcticament cap encavalcament.
Realitzem el mateix tipus de comparació per les variables FMps i ACps:
plot(num_norm[c(1,2)], col=ctg3clusters$cluster, main="Classificació k-means")
plot(num_norm[c(1,2)], col=as.factor(ctgp2$NSP), main="Classificació real")
legend("topright",
legend = c("normal", "sospitós", "patològic"),
col = 1:3, pch = 1)
plot(num_norm[c(1,2)], col=ctg2clusters$cluster, main="Classificació k-means")
plot(num_norm[c(1,2)], col=as.factor(ctgp2$NSP2), main="Classificació real")
legend("topright",
legend = c("normal", "sospitós"),
col = 1:2, pch = 1)
Tant en 2 com en 3 clústers, la divisió que proposa k-means difereix visiblement de les etiquetes diagnòstiques.
Fem una altra comparació, aquest cop només per 3 clústers, amb les variables ALTV i ACps:
plot(num_norm[c(1,7)], col=ctg3clusters$cluster, main="Classificació k-means")
plot(num_norm[c(1,7)], col=as.factor(ctgp2$NSP), main="Classificació real")
legend("topright",
legend = c("normal", "sospitós", "patològic"),
col = 1:3, pch = 1)
A la classificació real veiem com la categoria ‘normal’ es distribueix tot al llarg dels valors mitjans baixos de ACps, i desapareix per complet en els valors alts d’ALTV. En canvi, l’algorisme k-means situa dos grups diferents al llarg de la variable ACps, i també situa un tercer grup ocupant tots els valors de la variable ALTV que tenen 0 com a valors d’ACps.
Fem una última comparació entre les variables ALTV i FMps:
plot(num_norm[c(2,7)], col=ctg3clusters$cluster, main="Classificació k-means")
plot(num_norm[c(2,7)], col=as.factor(ctgp2$NSP), main="Classificació real")
legend("topright",
legend = c("normal", "sospitós", "patològic"),
col = 1:3, pch = 1)
Veiem un efecte semblant a la comparació anterior. A la classificació real, els valors baixos, propers a 0, de la variable FMps, mostren registres distribuïts a través de la variable ALTV de forma que podem identificar els tres grups diferents. El grup ‘normal’ té, o bé valor 0 de FMps, o bé valor 0 de ALTV. En canvi, a la classificació amb k-means, trobem un grup amb valors 0 de FMps, i els altres dos grups amb valors 0 d’ALTV. Es tracta doncs d’una classificació força diferent en el cas de la relació entre aquestes dues variables.
Ara generarem un segon model de classificació no supervisada, però aquest cop utilitzarem la fórmula de correlació de Pearson absoluta com a mètrica per mesurar la distància.
library(amap)
ctg3clusters_abpearson <- Kmeans(num_norm, 3, method = "abspearson")## Warning: did not converge in 10 iterationsComparem el nou model, amb 3 clústers, amb la classificació real, per les variables ALTV i ASTV:
plot(num_norm[c(6,7)], col=ctg3clusters_abpearson$cluster, main="Classificació k-means")
plot(num_norm[c(6,7)], col=as.factor(ctgp2$NSP), main="Classificació real")
Com en el model anterior, els 3 clústers que genera el model no s’acaben d’ajustar al tres grups diagnòstics.
Ara compararem, per les variables ALTV i FMps, el primer model que ja teníem, el model amb absolute Pearson, i la classificació real:
plot(num_norm[c(2,7)], col=ctg3clusters$cluster, main="Classificació k-means")
plot(num_norm[c(2,7)], col=ctg3clusters_abpearson$cluster, main="Classificació k-means/abspearson")
plot(num_norm[c(2,7)], col=as.factor(ctgp2$NSP), main="Classificació real")
legend("topright",
legend = c("normal", "sospitós", "patològic"),
col = 1:3, pch = 1)
Als gràfics veiem com, malgrat que els dos models k-means no són idèntics, tendeixen a dividir els clústers de forma similar, i cap dels dos s’aproxima a la classificació real.
Fem el mateix tipus de comparació per les variables ASTV i DPps, on veurem el mateix efecte:
plot(num_norm[c(5,6)], col=ctg3clusters_abpearson$cluster, main="Classificació k-means/abspearson")
plot(num_norm[c(5,6)], col=ctg3clusters$cluster, main="Classificació k-means")
plot(num_norm[c(5,6)], col=as.factor(ctgp2$NSP), main="Classificació real")
legend("topright",
legend = c("normal", "sospitós", "patològic"),
col = 1:3, pch = 1)
ctg2clusters_abpearson <- Kmeans(num_norm, 2, method = "abspearson")## Warning: did not converge in 10 iterations
Ara fem una última comparació per les mateixes variables ASTV i DPps, però aquest cop comparem per 2 clústers:
plot(num_norm[c(5,6)], col=ctg2clusters$cluster, main="Classificació k-means")
plot(num_norm[c(5,6)], col=ctg2clusters_abpearson$cluster, main="Classificació k-means/abspearson")
plot(num_norm[c(5,6)], col=as.factor(ctgp2$NSP2), main="Classificació real")
legend("topright",
legend = c("normal", "sospitós"),
col = 1:2, pch = 1)
El primer algorisme k-means talla horitzontalment el gràfic: divideix clarament els grups segons els valors de la variable ASTV. Amb la mètrica abspearson, l’algorisme classifica en un clúster els valors superiors de ASTV, i també aquells que tenen valors alts d’ambdues variables. En aquest sentit, la divisió de clústers s’assembla més a la classificació real, tot i que no classifica bé els valors intermedis.
if (!require('dbscan')) install.packages('dbscan'); library('dbscan')## Loading required package: dbscan##
## Attaching package: 'dbscan'## The following object is masked from 'package:fpc':
##
## dbscanFem una primera prova de l’algorisme DBSCAN amb un valor arbitrari d’eps = 1:
tryeps <- dbscan(num_norm, eps = 1)
tryeps## DBSCAN clustering for 2126 objects.
## Parameters: eps = 1, minPts = 5
## The clustering contains 1 cluster(s) and 0 noise points.
##
## 1
## 2126
##
## Available fields: cluster, eps, minPtsL’algorisme genera un sol clúster. Caldrà cercar un valor eps diferent.
Generem un gràfic de distància k que ens pugui donar pistes del valor òptim d’eps:
kNNdistplot(num_norm, k = 5)
abline(h=.25, col = "red", lty=2)
Amb la línia vermella hem indicat la zona on es podria trobar el valor d’eps que ens interessa, que seria entre 0.2 i 0.3.
Després de diferents intents entorn a aquests valors, trobem que amb el valor eps=0.29 obtenim 3 clústers i només 7 registres cauen fora de la classificació:
ctg_dbs <- dbscan(num_norm, eps = .29)
ctg_dbs## DBSCAN clustering for 2126 objects.
## Parameters: eps = 0.29, minPts = 5
## The clustering contains 3 cluster(s) and 7 noise points.
##
## 0 1 2 3
## 7 2088 23 8
##
## Available fields: cluster, eps, minPts
Ara utilitzarem l’algorisme OPTICS, i guardarem la variable que resulta d’extreure el dbscan per tal d’avaluar el model:
opti <- optics(num_norm)
opt <- extractDBSCAN(opti, eps_cl = .29)
opt## OPTICS ordering/clustering for 2126 objects.
## Parameters: minPts = 5, eps = 0.742812272170372, eps_cl = 0.29, xi = NA
## The clustering contains 3 cluster(s) and 7 noise points.
##
## 0 1 2 3
## 7 2088 8 23
##
## Available fields: order, reachdist, coredist, predecessor, minPts, eps,
## eps_cl, xi, cluster
Visualitzem els resultats de la classificació. Recordem que, també en la nostra classificació real, tenim un grup molt nombrós (diagnòstic ‘normal’) i dos grups molt menys nombrosos (diagnòstics ‘sospitós’ i ‘patològic’):
summary(as.factor(ctg1$NSP))## 1 2 3
## 1655 295 176Veiem que en els 3 clústers generats per l’algorisme optics aquesta diferència encara és més gran (només 8 i 23 registres en el 2n i 3r grup respectivament). Per aquest motiu, els gràfics no seran fàcils d’interpretar.
plot(opt)
Veiem una comparativa entre la classificació amb dbscan i la classificació real, per a totes les variables 2 a 2:
pairs(num_norm, col= opt$cluster)
pairs(num_norm, col= ctgp2$NSP)
Es fa difícil observar les diferències. Vegem doncs una comparació només per les variables ALTV i FMps:
plot(num_norm[c(2,7)], col=opt$cluster, main="Classificació DBSCAN")
plot(num_norm[c(2,7)], col=as.factor(ctgp2$NSP), main="Classificació real")
legend("topright",
legend = c("normal", "sospitós", "patològic"),
col = 1:3, pch = 1)
hullplot(num_norm[c(2,7)], opt)
Veiem en els gràfics que, comparant la classificació dbscan amb la classificació real, els resultats no semblen més gaire acurats respecte als models k-means que hem generat anteriorment.
Obtenim, amb l’índex de silhouette, una mètrica d’avaluació del model generat amb dbscan:
noise <- opt$cluster==0
clusters <- opt$cluster[!noise]
d <- dist(num_norm[!noise, 1:7])silh <- silhouette(clusters, d)
plot(silh, border=NA, col=sort(clusters), main="")
Si tenim en compte que els valors de l’índex de silhouette se situen entre -1 i 1, el valor que obtenim de 0.34 indica una qualitat d’agrupament acceptable.
En definitiva, els tres models comparats fins ara (kmeans, kmeans amb mètrica de Pearson, i dbscan) generen grups que s’emmirallen d’alguna manera amb la classificació real segons l’etiqueta diagnòstica NSP, sobretot per que fa als valors extrems dels parells de variables, i també en el sentit que agrupen 3 clústers de mides similars (1 grup nombrós i 2 grups minoritaris). Malgrat això, cap dels tres models encerta a separar de forma gaire acurada els registres amb combinacions de valors intermedis en les diferents variables.
Per la creació d’un model d’arbre de decisió, tenint en compte que l’aplicació i interpretació del model seran més fàcils i eficients si treballem amb variables categòriques, tenim la possibilitat de fer una nova selecció de variables. Tanmateix, per tal de seguir responent a l’objectiu inicial de recerca, ens interessa seguir incloent les variables que el sistema SisPorto mesura (recordem que no totes les variables de les quals disposem estan mesurades amb SisPorto). Així doncs, proposem fer algunes proves amb dos grups de variables. El primer inclourà la versió discretitzada de les variables mesurades amb SisPorto, a més de les variables diagnòstiques ‘AD’ (patró d’estrès) i ‘DE’ (patró de desacceleració), i el segon afegirà a aquesta mateixa selecció les variables diagnòstiques ‘A’ (son calmat), ‘B’ (son REM), ‘C’ (vetlla calmada), i ‘D’ (vetlla activa):
cat_ctg <- ctg1[c('AD', 'DE', 'AC_d', 'FM_d', 'UC_d', 'DL_d', 'DP_d', 'NSP2')]
cat_ctg_ext <- ctg1[c('A', 'B', 'C', 'D','AD', 'DE', 'AC_d', 'FM_d', 'UC_d', 'DL_d', 'DP_d', 'NSP2')]Creem conjunts train i test per a les dues mostres:
set.seed(666)
y <- cat_ctg[,8]
X <- cat_ctg[,1:7]set.seed(666)
y2 <- cat_ctg_ext[,12]
X2 <- cat_ctg_ext[,1:11]
split_prop <- 3
indexes = sample(1:nrow(cat_ctg), size=floor(((split_prop-1)/split_prop)*nrow(cat_ctg)))
trainX<-X[indexes,]
trainy<-y[indexes]
testX<-X[-indexes,]
testy<-y[-indexes]split_prop <- 3
indexes = sample(1:nrow(cat_ctg_ext), size=floor(((split_prop-1)/split_prop)*nrow(cat_ctg_ext)))
trainX2<-X2[indexes,]
trainy2<-y2[indexes]
testX2<-X2[-indexes,]
testy2<-y2[-indexes]
Comprovem que la divisió de conjunts sigui equilibrada:
sum(summary(trainy))## [1] 1417sum(summary(trainX$AD))## [1] 1417sum(summary(testy))## [1] 709sum(summary(testX$DE))## [1] 709sum(summary(trainy2))## [1] 1417sum(summary(trainX2$A))## [1] 1417sum(summary(testy2))## [1] 709sum(summary(testX2$B))## [1] 709
Mostrem les regles del primer arbre de decisió:
model1 <- C50::C5.0(trainX, trainy, rules=TRUE)
summary(model1)##
## Call:
## C5.0.default(x = trainX, y = trainy, rules = TRUE)
##
##
## C5.0 [Release 2.07 GPL Edition] Tue Jun 14 23:21:50 2022
## -------------------------------
##
## Class specified by attribute `outcome'
##
## Read 1417 cases (8 attributes) from undefined.data
##
## Rules:
##
## Rule 1: (713/18, lift 1.2)
## AC_d in {medium, high}
## DP_d in {low, medium}
## -> class Normal [0.973]
##
## Rule 2: (136/6, lift 1.2)
## DE = 1
## AC_d = low
## UC_d in {medium, high}
## DP_d in {low, medium}
## -> class Normal [0.949]
##
## Rule 3: (1015/98, lift 1.1)
## UC_d in {medium, high}
## DP_d = low
## -> class Normal [0.903]
##
## Rule 4: (8, lift 4.3)
## DE = 0
## AC_d = low
## DP_d = medium
## -> class Suspect [0.900]
##
## Rule 5: (68/8, lift 4.1)
## DP_d = high
## -> class Suspect [0.871]
##
## Rule 6: (186/59, lift 3.2)
## AC_d = low
## UC_d = low
## -> class Suspect [0.681]
##
## Default class: Normal
##
##
## Evaluation on training data (1417 cases):
##
## Rules
## ----------------
## No Errors
##
## 6 175(12.4%) <<
##
##
## (a) (b) <-classified as
## ---- ----
## 1050 67 (a): class Normal
## 108 192 (b): class Suspect
##
##
## Attribute usage:
##
## 87.09% DP_d
## 85.53% UC_d
## 73.61% AC_d
## 10.16% DE
##
##
## Time: 0.0 secsComparem amb l’arbre de decisió per la segona mostra:
model2 <- C50::C5.0(trainX2, trainy2, rules=TRUE)
summary(model2)##
## Call:
## C5.0.default(x = trainX2, y = trainy2, rules = TRUE)
##
##
## C5.0 [Release 2.07 GPL Edition] Tue Jun 14 23:21:50 2022
## -------------------------------
##
## Class specified by attribute `outcome'
##
## Read 1417 cases (12 attributes) from undefined.data
##
## Rules:
##
## Rule 1: (366, lift 1.3)
## B = 1
## -> class Normal [0.997]
##
## Rule 2: (275/2, lift 1.3)
## A = 1
## -> class Normal [0.989]
##
## Rule 3: (222/2, lift 1.3)
## AD = 1
## -> class Normal [0.987]
##
## Rule 4: (54, lift 1.3)
## D = 1
## -> class Normal [0.982]
##
## Rule 5: (40, lift 1.3)
## C = 1
## -> class Normal [0.976]
##
## Rule 6: (158/14, lift 1.2)
## DE = 1
## FM_d in {low, medium}
## -> class Normal [0.906]
##
## Rule 7: (300/1, lift 4.4)
## A = 0
## B = 0
## C = 0
## D = 0
## AD = 0
## DE = 0
## -> class Suspect [0.993]
##
## Rule 8: (2, lift 3.3)
## DE = 1
## FM_d = high
## -> class Suspect [0.750]
##
## Default class: Normal
##
##
## Evaluation on training data (1417 cases):
##
## Rules
## ----------------
## No Errors
##
## 8 19( 1.3%) <<
##
##
## (a) (b) <-classified as
## ---- ----
## 1097 1 (a): class Normal
## 18 301 (b): class Suspect
##
##
## Attribute usage:
##
## 47.00% B
## 40.58% A
## 36.84% AD
## 32.46% DE
## 24.98% D
## 23.99% C
## 11.29% FM_d
##
##
## Time: 0.0 secsEl primer arbre de decisió genera 5 regles, i obté una estimació d’error entorn un 11%. (* els percentatges d’errors i índexs de precisió que s’indiquen d’ara en endavant són aproximats, ja que varien lleugerament en cada execució de l’algorisme, i per tant també varien quan exportem el document de Rmd a Html).
Les variables més utilitzades són DP_d (desacceleracions perllongades), FM_d (moviments fetals), i AD (patró d’estrès), que és una variable diagnòstica.
Pel que fa a la segona mostra, amb 12 variables, l’arbre de decisió genera 7 regles, amb una estimació d’error de només entre un 1% i un 2%.
Mostrem gràficament els arbres dels dos models:
model1_plot <- C50::C5.0(trainX, trainy)
plot(model1_plot)
predicted_model1 <- predict( model1, testX, type="class" )
print(sprintf("La precisió de l'arbre de decisió es: %.4f %%",100*sum(predicted_model1 == testy) / length(predicted_model1)))## [1] "La precisió de l'arbre de decisió es: 84.4852 %"model2_plot <- C50::C5.0(trainX2, trainy2)
plot(model2_plot)
predicted_model2 <- predict( model2, testX2, type="class" )
print(sprintf("La precisió de l'arbre de decisió del segon model és: %.4f %%",100*sum(predicted_model2 == testy2) / length(predicted_model2)))## [1] "La precisió de l'arbre de decisió del segon model és: 98.5896 %"
D’aquests resultats, ens crida l’atenció la precisió del segon model. Sembla que és un model gairebé perfecte, amb una estimació d’error d’un 1.5% i una precisió de 98%. Això ens fa necessàriament sospitar de l’adequació del model. De fet, el primer model, amb un 86% de precisió, tampoc es queda curt.
La diferència entre els dos models té a veure amb la inclusió de variables diagnòstiques sumades a les variables mesurades per SisPorto. En el primer model en tenim dues, i en el segon model en tenim fins a 6. Una variable diagnòstica és, al capdavall, una variable classificatòria, i per tant té una funció similar a la variable NSP que estem utilitzant com a referència per supervisar el model. Pot ser doncs que això sigui un error de plantejament: si prenem com a exemple la variable AD (patró d’estrès), el sistema que classifica cada registre com a 0 o 1 en AD és el mateix sistema que classifica cada registre com a 1, 2 o 3 en la variable NSP, i és per tant poc sorprenent que la variable AD sigui predictora de NSP. Pot ser que el model sigui tan precís perquè reflecteix un efecte de circularitat: equivaldria a dir que un pacient està deprimit perquè té un trastorn depressiu, que és el mateix que dir que un pacient té un trastorn depressiu perquè està deprimit, quan el que ens interessa és establir que un pacient té un trastorn depressiu perquè té una desregulació del son, apatia, pensaments negatius, etc.
Així doncs, ens sembla que seria més correcte generar un model on prescindim de totes les variables classificatòries, i utilitzem només les variables mesurades amb SisPorto. Creiem que només així podem respondre al nostre objectiu inicial d’estudi. A més, això ens permetrà comparar els models no supervisats amb els models supervisats, ja que estarem treballant bàsicament amb les mateixes variables.
Hi ha dues variables mesurades per SisPorto, que hem utilitzat en els nostres models no supervisats, que encara no havíem discretitzat. Són ASTV (percentatge de temps amb variabilitat anormal de curta durada) i ALTV (percentatge de temps amb variabilitat anormal de llarga durada). Fem ara aquesta discretització:
library(arules)ctg_plus <- ctg1
ctg_plus$ASTV_d <- discretize(ctg_plus$ASTV, method = "cluster", labels = c("low", "medium", "high"))
ctg_plus$ALTV_d <- discretize(ctg_plus$ALTV, method = "cluster", labels = c("low", "medium", "high"))
summary(ctg_plus$ASTV_d)## low medium high
## 655 619 852summary(ctg_plus$ALTV_d)## low medium high
## 1704 275 147Creem ara una tercera mostra que inclou aquestes dues variables i descarta totes les variables diagnòstiques, excepte NSP2, que serà la nostra variable classificatòria:
cat_ctg3 <- ctg_plus[c('ASTV_d', 'ALTV_d','AC_d', 'FM_d', 'UC_d', 'DL_d', 'DP_d', 'NSP2')]
Fem els passos previs per generar l’arbre de decisió:
set.seed(666)
y3 <- cat_ctg3[,8]
X3 <- cat_ctg3[,1:7]split_prop <- 3
indexes = sample(1:nrow(cat_ctg3), size=floor(((split_prop-1)/split_prop)*nrow(cat_ctg3)))
trainX3<-X3[indexes,]
trainy3<-y3[indexes]
testX3<-X3[-indexes,]
testy3<-y3[-indexes]sum(summary(trainy3))## [1] 1417sum(summary(trainX3$ASTV_d))## [1] 1417sum(summary(testy3))## [1] 709sum(summary(testX3$ALTV_d))## [1] 709
Mostrem les regles d’aquest tercer arbre de decisió:
model3 <- C50::C5.0(trainX3, trainy3, rules=TRUE)
summary(model3)##
## Call:
## C5.0.default(x = trainX3, y = trainy3, rules = TRUE)
##
##
## C5.0 [Release 2.07 GPL Edition] Tue Jun 14 23:21:51 2022
## -------------------------------
##
## Class specified by attribute `outcome'
##
## Read 1417 cases (8 attributes) from undefined.data
##
## Rules:
##
## Rule 1: (854/49, lift 1.2)
## ASTV_d in {low, medium}
## -> class Normal [0.942]
##
## Rule 2: (1066/77, lift 1.2)
## ALTV_d = low
## DP_d in {low, medium}
## -> class Normal [0.927]
##
## Rule 3: (68/8, lift 4.1)
## DP_d = high
## -> class Suspect [0.871]
##
## Rule 4: (283/120, lift 2.7)
## ALTV_d in {medium, high}
## -> class Suspect [0.575]
##
## Default class: Normal
##
##
## Evaluation on training data (1417 cases):
##
## Rules
## ----------------
## No Errors
##
## 4 151(10.7%) <<
##
##
## (a) (b) <-classified as
## ---- ----
## 1066 51 (a): class Normal
## 100 200 (b): class Suspect
##
##
## Attribute usage:
##
## 95.20% ALTV_d
## 80.03% DP_d
## 60.27% ASTV_d
##
##
## Time: 0.0 secsL’arbre genera 4 regles de decisió, obté una estimació d’error d’un 10%, i utilitza tres variables: ALTV, DP i ASTV.
Mostrem gràficament l’arbre de decisió:
model3_plot <- C50::C5.0(trainX3, trainy3)
plot(model3_plot)
Veiem que aquest model obté una precisió del 87%:
predicted_model3 <- predict( model3, testX3, type="class" )
print(sprintf("La precisió de l'arbre de decisió es: %.4f %%",100*sum(predicted_model3 == testy3) / length(predicted_model3)))## [1] "La precisió de l'arbre de decisió es: 87.7292 %"
Ara mostrem en detall una matriu de confusió on podem comparar la predicció del model amb les dades del conjunt test:
if(!require(gmodels)){
install.packages('gmodels', repos='http://cran.us.r-project.org')
library(gmodels)
}## Loading required package: gmodels
CrossTable(testy3, predicted_model3,prop.chisq = FALSE, prop.c = FALSE, prop.r =FALSE,dnn = c('Reality', 'Prediction'))##
##
## Cell Contents
## |-------------------------|
## | N |
## | N / Table Total |
## |-------------------------|
##
##
## Total Observations in Table: 709
##
##
## | Prediction
## Reality | Normal | Suspect | Row Total |
## -------------|-----------|-----------|-----------|
## Normal | 517 | 21 | 538 |
## | 0.729 | 0.030 | |
## -------------|-----------|-----------|-----------|
## Suspect | 66 | 105 | 171 |
## | 0.093 | 0.148 | |
## -------------|-----------|-----------|-----------|
## Column Total | 583 | 126 | 709 |
## -------------|-----------|-----------|-----------|
##
##
Ara creem un segon model que inclogui adaptive boosting:
model3_ab <- C50::C5.0(trainX3, trainy3, trials = 100)
Veiem que amb aquesta modificació obtenim una petita millora (d’entorn un 1%) en la precisió del model:
predicted_model3_ab <- predict( model3_ab, testX3, type="class" )
print(sprintf("La precisió del model amb adaptive boosting: %.4f %%",100*sum(predicted_model3_ab == testy) / length(predicted_model3_ab)))## [1] "La precisió del model amb adaptive boosting: 88.7165 %"
Mostrem la matriu de confusió pel model amb adaptive boosting:
CrossTable(testy3, predicted_model3_ab, prop.chisq = FALSE, prop.c = FALSE, prop.r =FALSE,dnn = c('Reality', 'Prediction'))##
##
## Cell Contents
## |-------------------------|
## | N |
## | N / Table Total |
## |-------------------------|
##
##
## Total Observations in Table: 709
##
##
## | Prediction
## Reality | Normal | Suspect | Row Total |
## -------------|-----------|-----------|-----------|
## Normal | 511 | 27 | 538 |
## | 0.721 | 0.038 | |
## -------------|-----------|-----------|-----------|
## Suspect | 53 | 118 | 171 |
## | 0.075 | 0.166 | |
## -------------|-----------|-----------|-----------|
## Column Total | 564 | 145 | 709 |
## -------------|-----------|-----------|-----------|
##
## Si comparem les dues matrius de confusió, veiem que la millora de precisió es reflecteix prinicipalment en la detecció de veritables negatius (categoria ‘Suspect’). Respecta al primer model, el model amb adaptive boosting prediu més casos sospitosos que realment ho són. Al mateix temps, el nombre de casos classificats erròniament com a normals es redueix amb adaptive boosting. Així doncs, malgrat que l’augment de precisió és només d’entorn els 2 punts percentuals, el segon model respon de forma més efectiva als nostres objectius.
Ara generem un model amb l’algorisme Random Forest que classifica els registres segons les mateixes variables que hem utilitzat fins ara. Dividim els conjunts d’entrenament i test:
set.seed(123)
samp <- sample(nrow(cat_ctg3), 0.6669 * nrow(cat_ctg3))
train <- cat_ctg3[samp, ]
test <- cat_ctg3[-samp, ]
I apliquem el model:
library(randomForest)## randomForest 4.7-1## Type rfNews() to see new features/changes/bug fixes.##
## Attaching package: 'randomForest'## The following object is masked from 'package:ggplot2':
##
## marginmodel_rf <- randomForest(NSP2~., data = train)
model_rf##
## Call:
## randomForest(formula = NSP2 ~ ., data = train)
## Type of random forest: classification
## Number of trees: 500
## No. of variables tried at each split: 2
##
## OOB estimate of error rate: 9.74%
## Confusion matrix:
## Normal Suspect class.error
## Normal 1092 23 0.0206278
## Suspect 115 187 0.3807947Obtenim una estimació d’error del 11%, havent utilitzat 2 variables a cada bifurcació del model.
Amb aquest model obtenim una precisió del 88%, un nivell similar als models d’arbre de decisió anteriors (lleugerament per sota del model amb adaptive boosting):
prediction_rf <- predict(model_rf, newdata = test)
sum(prediction_rf==test$NSP2) / nrow(test)## [1] 0.8984485A la representació gràfica del model veiem com les variables que es consideren com a més importants en la predicció del model són ASTV, ALTV, AC i DP (el nombre ’MeanDecreaseGini indica quin percentatge decauria la precisió de la predicció si alguna d’aquestes variables no es tingués en compte):
varImpPlot(model_rf)
Tot seguit, mostrem les matrius de confusió dels models d’arbre decisió sense i amb boosting, i del model Random Forest:
library(caret)## Loading required package: lattice
confusionMatrix(data= predicted_model3, reference = testy3)## Confusion Matrix and Statistics
##
## Reference
## Prediction Normal Suspect
## Normal 517 66
## Suspect 21 105
##
## Accuracy : 0.8773
## 95% CI : (0.8509, 0.9005)
## No Information Rate : 0.7588
## P-Value [Acc > NIR] : 1.772e-15
##
## Kappa : 0.6317
##
## Mcnemar's Test P-Value : 2.390e-06
##
## Sensitivity : 0.9610
## Specificity : 0.6140
## Pos Pred Value : 0.8868
## Neg Pred Value : 0.8333
## Prevalence : 0.7588
## Detection Rate : 0.7292
## Detection Prevalence : 0.8223
## Balanced Accuracy : 0.7875
##
## 'Positive' Class : Normal
##
confusionMatrix(data=predicted_model3_ab, reference = testy3)## Confusion Matrix and Statistics
##
## Reference
## Prediction Normal Suspect
## Normal 511 53
## Suspect 27 118
##
## Accuracy : 0.8872
## 95% CI : (0.8615, 0.9095)
## No Information Rate : 0.7588
## P-Value [Acc > NIR] : < 2.2e-16
##
## Kappa : 0.6749
##
## Mcnemar's Test P-Value : 0.005189
##
## Sensitivity : 0.9498
## Specificity : 0.6901
## Pos Pred Value : 0.9060
## Neg Pred Value : 0.8138
## Prevalence : 0.7588
## Detection Rate : 0.7207
## Detection Prevalence : 0.7955
## Balanced Accuracy : 0.8199
##
## 'Positive' Class : Normal
##
confusionMatrix(data=prediction_rf, reference = test$NSP2)## Confusion Matrix and Statistics
##
## Reference
## Prediction Normal Suspect
## Normal 534 66
## Suspect 6 103
##
## Accuracy : 0.8984
## 95% CI : (0.8738, 0.9197)
## No Information Rate : 0.7616
## P-Value [Acc > NIR] : < 2.2e-16
##
## Kappa : 0.6815
##
## Mcnemar's Test P-Value : 3.57e-12
##
## Sensitivity : 0.9889
## Specificity : 0.6095
## Pos Pred Value : 0.8900
## Neg Pred Value : 0.9450
## Prevalence : 0.7616
## Detection Rate : 0.7532
## Detection Prevalence : 0.8463
## Balanced Accuracy : 0.7992
##
## 'Positive' Class : Normal
## Ja havíem comparat els dos primers models, i vèiem que el model amb adaptive boosting obté major precisió, però sobretot obté una especificitat molt més gran (classifica correctament com a casos sospitosos els que realment ho són). El model Random Forest no millora els resultats del model amb adaptive boosting, tot i que la seva sensibilitat és molt alta.
El nostre estudi tenia com a objectiu plantejar models amb una alta capacitat de classificació. Concretament, ens interessava optimitzar l’especificitat dels possibles models de selecció, ja que imaginàvem el cas en què en l’anàlisi d’algunes cardiotocografies, el sistema SisPorto havia diagnosticat com a normal alguns fetus amb patrons anòmals de ritme cardíac.
El nostre joc de dades tenia fins a 40 variables de diferents tipus, incloent les mesures del sistema SisPorto, però també altres variables amb càlculs derivats dels registres, i variables categòriques diagnòstiques.
Per tal de centrar-nos en l’objectiu d’anàlisi, i també per poder comparar els resultats dels diferents models, ens hem cenyit a les variables mesurades pel sistema SisPorto tant pels models no supervisats com pels models supervisats (hem discretitzat variables per adaptar-les als models supervisats).
Hem constatat que els models no supervisats no mostren, en el cas del nostre joc de dades, una gran capacitat de classificar correctament segons els 3 grups diagnòstics (variable NSP: ‘normal’, ‘sospitós’ i ‘patològic’). Val a dir però que els models no supervisats k-means sí que identifica un nombre òptim de clústers (2 o 3) que té relació amb aquesta mateixa variable classificatòria que nosaltres coneixem d’entrada. Quan, per avaluar els models no supervisats que hem generat (k-means, k-means amb mètrica absolute Pearson i DBSCAN), fem comparacions entre l’agrupació de l’algorisme i les dades reals per 2 variables, veiem repetidament que els models tenen dificultats per discriminar correctament els registres que mostren valors intermedis de les variables. Una altra dificultat, que es fa visible especialment en l’aplicació del model DBSCAN, és que els nostres 3 grups diagnòstics tenen mides molt desiguals: els casos patològics representen una part petita dels registres, i la majoria de casos són diagnosticats com a normals. El model DBSCAN exagera aquesta desigualtat, i classifica com a patològics una quantitat excessivament petita de registres. En definitiva, un model no supervisat, almenys entre els que hem estudiat, no seria una bona opció per detectar casos patològics de forma automàtica.
Pel que fa als models supervisats, hem trobat que, en general, aconseguíem uns índexs de precisió molt alts. La quantitat de registres i variables semblava formar una bona base per l’entrenament dels algorismes, tant pel que fa als arbres de decisió com l’algorisme Random Forest. En ambdós casos, no semblava tampoc que hi hagués un problema de sobreentrenament, ja que els algorismes seguien mostrant bons resultats sobre els conjunts test. Hem constatat que, entre els 3 models supervisats que hem generat, el model d’arbre de decisió amb adaptive boosting mostra el nivell d’especificitat més alt, i classifica erròniament casos patològics amb menor freqüència comparat amb els altres dos models. Tanmateix, si imaginem utilitzar aquest algorisme per detectar de forma automàtica els casos patològics, trobem que encara hi ha massa casos que es classifiquen com a normals que en realitat són patològics. I aquí rau el perill d’utilitzar aquests models en el nostre cas concret. Tot i que hem pogut generar models supervisats amb bons índexs de precisió, cal tenir present la rellevància de l’error que suposa diagnosticar com a normal una cardiotocografia d’un fetus que té un patró patològic de ritme cardíac. Això no significa necessàriament que els nostre models, o d’altres més òptims que es puguin generar, no siguin útils, sinó que cal seguir utilitzant-los com a suport o contrast a l’avaluació d’un expert que pugui tenir en compte altres variables més enllà de les mesurades per sistema automàtic SisPorto.